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We express the transport coefficients appearing in the second-order evolution eqnations for bulk 
viscous pressure and shear stress tensor using Bose-Einstein, Boltzmann, and Fermi-Dirac statistics 
for the equilibrium distribution function and Grad’s 14-moment approximation as well as the method 
of Chapman-Enskog expansion for the non-equilibrium part. Specializing to the case of transversally 
homogeneous and boost-invariant longitudinal expansion of the viscous medium, we compare the 
results obtained using the above methods with those obtained from the exact solution of the massive 
0-|-ld relativistic Boltzmann equation in the relaxation-time approximation. We show that compared 
to the 14-moment approximation, the hydrodynamic transport coefficients obtained by employing 
the Chapman-Enskog method leads to better agreement with the exact solution of the relativistic 
Boltzmann equation. 
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I. INTRODUCTION 

Relativistic viscous hydrodynamics has been applied 
quite successfully to study and understand various collec¬ 
tive phenomena observed in the evolution of the strongly 
interacting QCD matter, with very high temperature and 
density, created in relativistic heavy-ion collisions; see 
Ref. [1] for a recent review. The derivation of hydrody¬ 
namic equations is essentially a coarse graining procedure 
whereby one obtains an effective theory describing the 
long-wavelength low-frequency limit of the microscopic 
dynamics of a system [1, 2]. Relativistic hydrodynamics 
is formulated as an order-by-order expansion in powers 
of space-time gradients where ideal hydrodynamics is of 
zeroth order [1]. The viscous effects arising in the first- 
order theory, also known as the relativistic Navier-Stokes 
theory [2, 3], results in acausal signal propagation and 
numerical instability. While causality is restored in the 
second-order Israel-Stewart (IS) theory [4], stability may 
not be guaranteed [5] . Consistent formulation of a causal 
theory of relativistic viscous hydrodynamics and accurate 
determination of the associated transport coefficients is 
currently an active research topic [6-33]. 

The second-order IS equations can be derived in sev¬ 
eral ways [5]. For example, in the derivations based on 
the second law of thermodynamics, the hydrodynamic 
transport coefficients related to the relaxation times for 
bulk and shear viscous evolution remain undetermined. 
While these transport coefficients can be obtained in the 
derivations based on kinetic theory [4, 6] , the form of non¬ 
equilibrium phase-space distribution function, f(x,p), 
has to be specified. Two most extensively used meth¬ 
ods to determine f{x,p) for a system which is close 
to local thermodynamic equilibrium are (1) the Grad’s 
14-moment approximation [34] and (2) the Chapman- 
Enskog method [35]. Note that while Grad’s 14-moment 


approximation has been widely employed in the formu¬ 
lation of a causal theory of relativistic dissipative hydro¬ 
dynamics [4-14], the Chapman-Enskog method remains 
less explored [15-19]. On the other hand, the Chapman- 
Enskog formalism has been often used to extract vari¬ 
ous transport coefficients of hot hadronic matter [20-24] 
Although in both methods the distribution function is 
expanded around its equilibrium value fo{x,p), it has 
been demonstrated that the Chapman-Enskog method in 
the relaxation-time approximation (RTA) leads to better 
agreement with both microscopic Boltzmann simulations 
as well as exact solutions of the relativistic RTA Boltz¬ 
mann equation [16-19]. This may be attributed to the 
fact that a fixed-order moment expansion, as required in 
Grad’s approximation, is not necessary in the Chapman- 
Enskog method. 

Much of the research on the application of viscous hy¬ 
drodynamics in relativistic heavy-ion collisions is devoted 
to the extraction of the shear viscosity to entropy den¬ 
sity ratio, 77 /s, from the analysis of the anisotropic flow 
data [36-38]. Indeed the estimated 77 /s has been found to 
be close to the conjectured lower bound 77/sjKSS = l/ 47 r 
[39, 40]. On the other hand, a self-consistent and sys¬ 
tematic study of the effect of the bulk viscosity in nu¬ 
merical simulations of high-energy heavy-ion collisions is 
relatively lacking. This may be attributed to the fact that 
the hot QCD matter is assumed to be nearly conformal 
and therefore the bulk viscosity is estimated to be much 
smaller compared to the shear viscosity. However, in re¬ 
ality, QCD is a non-conformal field theory and therefore 
bulk viscous corrections to the energy-momentum tensor 
should not be ignored in order to correctly understand 
the dynamics of the QCD system. Moreover, for the 
range of temperature explored in heavy-ion collision ex¬ 
periments at Relativistic Heavy-Ion Collider (RHIC) and 
Large Hadron Collider (LHC), the magnitude and tern- 
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perature dependence of the bulk viscosity could be large 
enough to influence the space-time evolution of the hot 
QCD matter [41-47]. 

It is important to note that the second-order transport 
coefficients, appearing in the evolution equation for the 
bulk viscous pressure, are less understood compared to 
those of the shear stress tensor. While the relaxation 
time for the bulk viscous evolution can be obtained by 
using the second law of thermodynamics in a kinetic 
theory set up [11, 12], this method fails to account for 
the important coupling of the bulk viscous pressure with 
the shear stress tensor [13, 14, 18]. On the other hand, 
for finite masses and classical Boltzmann distribution, 
the second-order transport coefficients corresponding to 
bulk viscous pressure and shear stress tensor have been 
obtained by employing the Grad’s 14-moment approxi¬ 
mation [13, 14] as well as the Chapman-Enskog method 
[18], within a purely kinetic theory framework. However, 
these transport coefficients still remain to be determined 
for quantum statistics, i.e., for Bose-Einstein and Fermi- 
Dirac distribution. 

In this paper, we express the transport coefficients ap¬ 
pearing in the second-order viscous evolution equations 
with non-vanishing masses for Bose-Einstein, Boltzmann 
and Fermi-Dirac distribution. We obtain these trans¬ 
port coefficients using the Grad’s 14-moment approxi¬ 
mation as well as the method of Ghapman-Enskog ex¬ 
pansion. In addition, in the case of one-dimensional 
scaling expansion of the viscous medium, we compare 
the results obtained using the above methods with those 
obtained from the exact solution of massive 0-|-ld rel¬ 
ativistic Boltzmann equation in the relaxation-time ap¬ 
proximation [33]. We demonstrate that the results ob¬ 
tained using the Ghapman-Enskog method are in bet¬ 
ter agreement with the exact solution of the RTA Boltz¬ 
mann equation than those obtained using the Grad’s 14- 
moment approximation. 


II. RELATIVISTIC HYDRODYNAMICS 

In the absence of any conserved charges, i.e., for van¬ 
ishing chemical potential, the hydrodynamic evolution of 
a system is governed by the local conservation of energy 
and momentum: = 0. The energy-momentum 

tensor, , which characterizes the macroscopic state of 
a system, can be expressed in terms of the single-particle 
phase-space distribution function, f{x,p), and tensor de¬ 
composed into hydrodynamic degrees of freedom [48], 

= JdPp^p''f{x,p) = eu^^u'' -{P + -t 

( 1 ) 

Here dP = gd^p/[{2Tr)^p^] is the invariant momentum- 
space integration measure, g being the degeneracy fac¬ 
tor, and p^ is the particle four-momentum. In the tensor 
decomposition, e, P, H, and are the energy den¬ 
sity, the thermodynamic pressure, the bulk viscous pres¬ 


sure, and the shear stress tensor, respectively. The pro¬ 
jection operator = g^'' — is constructed such 

that it is orthogonal to the hydrodynamic four-velocity 
. The metric tensor is assumed to be Minkowskian, 
g^^'^ = diag(-|-l, —1, —1, —1), and is defined in the Lan¬ 
dau frame: 

The energy-momentum conservation equation, 
= 0, when projected along and orthogonal to 
gives the evolution equations for e and u^, 

e {e P n)0 — (J= 0, 
(e-bP-bn)M“-V“(P-bn)-bA“9^7r'"'' = 0. (2) 

Here we have used the usual notation 9 = for the 

expansion scalar, A = u^d^A for the co-moving deriva¬ 
tive, V“ = A'^“9^ for the space-like derivative, and 
+ V‘^u'")/2 - (0/3)A'"'^ for the velocity 
stress tensor. The energy density and the thermody¬ 
namic pressure can be written in terms of the equilibrium 
phase-space distribution function /q as 


eo = u^Ui, J dPp^p'' fo, 

(3) 

Po = -^A^.JdPp^p’'fo, 

(4) 


where the equilibrium distribution function for vanishing 
chemical potential is given by 


exp (/3m ■ p) + a 

Here u ■ p = u^p^ and a = —1,0,1 for Bose-Einstein, 
Boltzmann, and Fermi-Dirac gas, respectively. The in¬ 
verse temperature, /3 = l/T, is determined by the match¬ 
ing condition e = cq. 

For a system close to local thermodynamic equilib¬ 
rium, the non-equilibrium phase-space distribution func¬ 
tion can be written as / = fo + Sf, where 6f /. Using 
Eq. (1), the bulk viscous pressure, H, and the shear stress 
tensor, can be expressed in terms of Sf as [48] 

n = -^A^pJdPpyf^6f, ( 6 ) 

7r'^^ = A^^y'dPpV<5/, (7) 

where A^^ = (A([A^ + A^^A^)/2 - (1/3)A^"A„;3 is a 
traceless symmetric projection operator orthogonal to 
and . In the following, we employ the expressions for 
Sf obtained using the Grad’s 14-moment approximation 
and the Chapman-Enskog like iterative solution of the 
relativistic Boltzmann equation to obtain expressions for 
the quantum transport coefficients associated with vis¬ 
cous evolution. 


III. VISCOUS CORRECTIONS TO THE 
DISTRIBUTION FUNCTION 

Precise determination of the form of the non¬ 
equilibrium single particle phase-space distribution func- 
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tion is one of the central problems in statistical physics. 
For a system close to local thermodynamic equilibrium, 
the problem reduces to determining the form of the cor¬ 
rection to the equilibrium distribution function. Within 
the framework of relativistic hydrodynamics, the viscous 
corrections to the equilibrium distribution function can 
be obtained from two different methods: (1) the moment 
method and (2) the Chapman-Enskog method. The mo¬ 
ment method, more popularly known as the Grad’s 14- 
moment ansatz, is based on a Taylor-like expansion of 
the non-equilibrium distribution in powers of momenta. 
On the other hand, the Chapman-Enskog method relies 
on the solution of the Boltzmann equation. 

Ignoring dissipation due to particle diffusion, the 
Grad’s 14-moment approximation leads to 


Sfo 


{Eq + B^rn^ E Do{u ■ p) - 4Bo(u ■ pY] 11 


and Cg = dPjde is the speed of sound squared which can 
be expressed as 


c 


2 

S 


e + P 

p4o’ ■ 


(17) 


IV. VISCOUS EVOLUTION EQUATIONS 

The second-order evolution equations for 11 and 
can be derived by considering the co-moving derivative 
of Eqs. (6) and (7), 

n = -^A^pJdPpy^Sf, (18) 

=A^J'^JdPp^p^Sf. (19) 


B2P°‘p'^na/3 


fofo, 


( 8 ) 


where /o = 1 — a/o- The coefficients Eq, Bq, Dq and B 2 
are known in terms of m, T and u-p and can be expressed 
as 


B. = 


( 0 )’ 


3So 

Eq 


9 T 

t(0) t(0) _ t(0) j (0) 

_ I ’^31 ’^20 ’^41 ’^10 _ ^ 


t (0) t (0) 

^30 "^10 


_ t(0) t(0) 

’^20 "^20 


t(0) t(0) _ 7 ( 0 ) t(0) 

^2 _ ^ "^31 *^30 ’^41 '^20 


j: 


(0) r(0) 


30 


Ti 

•^10 "^20 '^20 


Bo=- 


(0) ,^(0) 

20 

1 


= Ci, 


-I- qC t -I- ^ 

oL^lJ2i ~r 0^2^31 — '^'^41 “r UJ 42 


(9) 


( 10 ) 


( 11 ) 


( 12 ) 


(r) 

where the thermodynamic functions Jnq are defined as 


^ (2gTl)!! 

(13) 

In the above equation, the indices n — r and q represents 
the number of times p^ and appear in the integra¬ 
tion, respectively. 

An analogous expression for 6f can also be obtained 
using an iterative Chapman-Enskog like solution of the 
relativistic Boltzmann equation in the relaxation-time 
approximation [18, 49]. In absence of dissipation due 
to particle diffusion, the Chapman-Enskog method leads 
to [18] 


SfcE = 


J3_ 

u ■ p 


1 

3/3n 

1 


{m^- (l-3cf)(u-p)^}n 

fofo, (14) 


+ 2^P>X. 


where 

l3u = lp.-{e + P)cl (15) 

I3^ = (3jg\ (16) 


Further, Sf appearing in the above equations can be sim¬ 
plified by rewriting the relativistic Boltzmann equation, 

P^d^f = C[f], (20) 

in the form [9] 

Sf = -fo-—p^V^f+ — C[f]. (21) 

u-p u-p 

In the following, we consider the relaxation-time approx¬ 
imation for the collision term in the Boltzmann equation, 

C[f] = -{u-p)^, (22) 

where Teq is the relaxation time. In general, in order 
for the collision term in the above equation to respect 
the conservation of particle four-current and energy- 
momentum tensor, Teq must be independent of momenta 
and u^ has to be defined in the Landau frame [50]. 

Substituting Sf from Eq. (21) into Eqs. (18) and (19) 
along with the form of Sf given in Eqs. (8) and (14), and 
after performing the integrations, we obtain 

n =- /3nd — (5nnII6< -I- (23) 

T ' 

' eq 

- (24) 

where — \/''u^) is the vorticity tensor. 

Note that the above form of the evolution equations for 
the bulk viscous pressure and the shear stress tensor 
is exactly same for both Grad’s 14-moment approxima¬ 
tion {Sfc) and the Chapman-Enskog expansion (Sfcs)- 
Moreover, the expressions for the first order transport co¬ 
efficients, /3n and /3,r, are also identical for these two cases 
and are given in Eqs. (15) and (16). However, the second 
order transport coefficients appearing in the above equa¬ 
tions are different for the Grad’s 14-moment method and 
the Chapman-Enskog method. 














4 


The transport coefficients in the case of Grad’s 14- 
moment approximation are calculated to be 


- 1 _ c2 _ 

"nn — ^ G n '2 ) 


A(G) _ i _ 

— 3 G + 3 72 I 

J(G) _ 4 

^TTTT 3 ^ 3 ’ 


r(G) = 

TTTT 


10 4m^ 


7 


7 


-72 


( 2 ) 


(G) _ 6 2m4 (0) 

- 5 15 72 , 


(25) 

(26) 

(27) 

(28) 
(29) 


where 

„(o) 


7^' = (Eo + Bom^)J% + DoJ^’lo - 4i?o (30) 


(0) 


r(0) 


„( 2 ) 


72 = 


( 0 ) 

22 


7(0)' 

^42 


(31) 


On the other hand, these transport coefficients in the 
case of Chapman-Enskog method are obtained as 


a(CE) 

"nn 

5 2 

= -gX-G, 

(32) 

7CE) 

— ^ (7 -1-2 7^^^) — 

“ 3/3^ V + 4^42 j C„ 

(33) 

J(CE) 

^TTTT 

_ 5 , 7/3 (3) 

3 3/3^ ’ 

(34) 

^(CE) 

' TTTT 

— 2 -1- — 

-2+/3^763, 

(35) 

7CE) 

2 

= -3%, 

(36) 


where 


X = 


/3n 


(1 - ici)[jg'> + Jg^) - + 4g) 


(37) 

(r) 

The integral functions Jnq appearing in the expres¬ 
sions for the transport coefficients satisfy the relations 


= 


2 j(r) 


m J' 


-J 


(r) 


"9 (2g-I- 1) L "-2.9-1 '^",9-1 


/(O) = 1 




(0) 


(38) 

(39) 


where. 


/(") = _ 

“ (2g+l)!!. 


dP {u ■ p) 


n—2q—r 


{A^^pVrfo. (40) 


Here we readily identify = e and gg = —P- Using 
Eqs. (38) and (39), we obtain the identities 


(0) _ rr4_ (0) _ 1 70) 
•^42 “ 5 '^21 5'^41 J 

jgi=-\^^+n 


(41) 

(42) 


To compute all the transport coefficients, we also need to 
determine the integrals J 2 o\ jgi ^ '^ 3 o^ '4i'20J 

J-L and Jg4 In the follow¬ 

ing, we obtain expressions for these quantities in terms 
of modified Bessel functions of the second kind. 


V. TRANSPORT COEFFICIENTS 


Let us first simplify our equilibrium distribution func¬ 
tion. Using the result of summation of a infinite geomet¬ 
ric progression, 

. CXD 

- = \ — x + x^ — x^--- = (—l)*a;*, for |x| < 1, 

' + »■ ES 

(43) 

we obtain. 


g-^("'p) 

1 -I- ae-^(“-P) 


1^1 


(44) 


Similarly, using the result obtained after differentiating 
Eq. (43), we obtain 


/o/o 


g-/3(«.p) 

(1 -f ae-^(“'P))2 


(45) 

/ = 1 


Using the above results for /o and /o/o, the thermody- 

(r^ fr) 

namic integrals Inq and Jnq can be written as. 


r(") 

^nq 


/(") 
^ nq 


g'j'n+2 — r ^n-\-2 — r 

27r2(2g-h 1)!! 

X (cosh6l)”-2«-’’ 

g'j'n-\-2 — r ^n-\-2 — r 


27r2(2g-h 1)!! 

X (cosh6l)”-2«-’’ 


dO (46) 

i=i -lo 

(sinh0)^'^-'-^ exp(—Zzcosh0), 

pGO 

{-iyY.^{-ay~4 dO (47) 
i=i -lo 

(sinh0)^'^-'-^ exp(—Zzcosh0), 


where z = Pm = m/T. 

(r) fr) 

The integral coefficients Inq and Jnq can be expressed 
in terms of modified Bessel functions of the second kind. 
The integral representation of the relevant Bessel func¬ 
tion is given by 


Kn{z) = / d6*cosh(n0) exp(—zcosh0). (48) 

Jo 

The pressure, energy density and Jgg^, required to calcu¬ 
late the velocity of sound, can be expressed in terms of 
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Bessel functions using the above equation 
“ 1 


P = 


e = 


1 = 1 

1 


27r2 I 

1^1 


J2 (-ay-^K,{lz)-P, 


(49) 

(50) 


4o = F5(Z^) + Ksilz) - 2K,{lz) 


1^1 


(51) 


Here the function Ki^n{lz) is defined by the integral 

poo in 

K,^n{lz)=j^ (coshg)" exp(-^zcoshd). (54) 

Note that the subscript i in the above function is not an 
index and just serves to distinguish it from the Bessel 
functions. Ki^n satisfy the following recurrence relation: 

■^Ki^nilz) = -lKi^n-i{lz), (55) 

which can be expressed in the integral form: 


(r) 

The integral functions Jnq appearing in the transport 
coefficients obtained using the Grad’s 14-moment method 
can be expressed as 




(0)_ gT^z^f.1 

■^21 - 271^ ^ I 

1^1 

( 0 ) _ gT z -i- 
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7 ( 0 ) _ 
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1=1 
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1927r2 ■ 






Kg - 2Ka -K 2 + 2Kn 
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Z=1 
p4 4 00 


Z=1 


2^2 


7 ( 0 ) _ 9 T Z 
■^00 - 4^2 


E^(-F-' 


1=1 


K.-Kn 


(52) 


where the (^ 2 :)-dependence of and Ki ^ is implic¬ 
itly understood. The integral functions appearing in 
the transport coefficients obtained using the Chapman- 
Enskog method are 


j(3) __ 

"'63 - 


33607r2. 


K 5 -IIK^P 58Ki 


1=1 


- 641^,4 + 16K,,3 
,5 

“ 480772 


4 y = £^^Z(-a)'-i K5- 7X3+ 221 ^ 1 - 16 i^,.i 


Z = 1 


7(3) _ 
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gTh 


I 2 O 772 

(2) _ gT^z' 

t/QI - 


^/(-a)'-i F3- 9X1+ 12 i^,,i- 4 i ^,,3 


Z=1 
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247r2 


7^3-5Ki+4X,,i . (53) 


1=1 


K,^n{lz) = if*.„(0) - I [ K,^n-i{lz')dz'. (56) 

Jo 

We observe that by matching Kifi(lz) = Ko{lz), the 
above recursion relation can be used to evaluate Ki^n(lz) 
up to any n. 

( 7 *) 

Armed with the above expressions for Jnq in terms 
of series summation of the Bessel function, it is instruc¬ 
tive to calculate the ratio between the coefficient of bulk 
viscosity and shear viscosity. In the relaxation-time ap¬ 
proximation, this ratio is given by (/g = PuIP-k [13, 18]. 
In the small-z approximation, using the series expansion 
of the Bessel functions in powers of 2 , we obtain 

^ = +0{z^), (57) 


where T = 75 for Boltzmann distribution and T = 48 
for Fermi-Dirac distribution. However, in the case of 
Bose-Einstein distribution (remember a = —1), we get 
r = —15 -I- 36/(1 -I- a), and therefore it diverges as 
~ 1/(1 -I- a) up to the leading-order. It is interesting 
to note that the above expression is similar to the well 
known relation, (/g = 15(1/3 — c^)^, derived by Wein¬ 
berg [51]. The difference in the proportionality constant 
may be attributed to the fact that here we have consid¬ 
ered a single component system whereas, motivated by 
applications to cosmology, Weinberg considered a system 
composed of mixture of radiation and matter. 

In order to calculate (^/g for the quark-gluon plasma 
(QGP), we need to provide appropriate degeneracy fac- 

(r) ('r') 

tors for the integral coefficients Ikq and Jnq . For QGP, 
these integral coefficients will transform as 



g=gg,a=i 

g=gq,a=i 


■ - - 1 

^nq \9—9g,Q‘— — l'> 

- 

'■'ng \g=gg,a=-lj 


(58) 


where gq and gg are the degeneracy factors corresponding 
to quarks and gluons respectively. These factors are given 
as 


gg= NgX Nqg X Nc X Nf = 12Nf, 
gg=NgX {Nl - 1) = 16, (59) 

where Nf is the number of flavours, Nc = 3 is the num¬ 
ber of colours, Ns = 2 corresponds to spin degrees and 
Nqq = 2 denotes quark and anti-quark. 
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FIG. 1: (Color online) mjT dependence of the ratio C,lr\ 
scaled by (1/3 — for Boltzmann (red solid line), Bose- 
Einstein (blue dashed line) and Fermi-Dirac (green dotted 
line) statistics. We also show results for two-flavor (brown 
dashed-dotted line) and three-flavor QGP (purple dashed- 
dotted-dotted line). Thin horizontal black dashed lines cor¬ 
responds to constant values 75 and 48. 


In Fig. 1, we show m/T dependence of the ratio C/t? 
scaled by (1/3—c^)^ for various cases. We observe that in 
accordance with the predictions from small-z expansion, 
Eq. (57), the curve for Boltzmann statistics (red solid 
line) and Fermi-Dirac statistics (green dotted line) satu¬ 
rates at 75 and 48 (thin horizontal black dashed lines), 
respectively, at small values of z. We also see that at 
very small- 2 , Bose-Einstein statistics (blue dashed line) 
results in very large values of the scaled viscosity ratio 
(C/^)/(l/3 —Cs)^ indicating divergence. We see that this 
ratio for the QGP is dominated by the Bose-Einstein 
statistics for two-flavor (brown dashed-dotted line) as 
well as three-flavor QGP (purple dashed-dotted-dotted 
line). This however does not imply that the bulk viscos¬ 
ity of QGP is very large because a realistic calculation 
with lattice QGD equation of state suggest that z > 0.6 
even for temperatures as high as central RHIG and LHG 
energies [49, 52, 53]. Moreover, for large-z we see that 
all three statistics lead to similar results for the scaled 
viscosity ratio which suggests that there is a relatively 
small effect coming from composition of the fluid. 


VI. BOOST-INVARIANT 0-blD CASE 


In this section, we consider evolution in the case of 
transversely homogeneous and purely-longitudinal boost- 
invariant expansion [54] . For such an expansion all 
scalar functions of space and tim e depen d only on the 
longitudinal proper time r = y/t'^ — z'^. Working in 
the Milne coordinate system, {T,x,y,ri), the hydrody¬ 
namic four-velocity is given by = (1,0,0,0). The 
energy-momentum conservation equation together with 



T [fm/c] 


FIG. 2: (Color online) Time evolution of the the bulk viscous 
pressure times r (top) and the pressure anisotropy Pl/Pt 
scaled by that obtained using exact solution of the RTA Boltz¬ 
mann equation (bottom) for Boltzmann statistics. The three 
curves in both panels correspond to three different calcula¬ 
tions: the exact solution of the RTA Boltzmann equation [33] 
(red solid line), second-order viscous hydrodynamics obtained 
using the Chapman-Enskog method (blue dashed line) and us¬ 
ing the Grad’s 14-moment approximation (green dotted line). 
For both panels we use Tq = 300 MeV at tq = 0.5 fm/c, 
m = 300 MeV, and Teq = — rn = 0.25 fm/c. The initial 

spheroidal anisotropy in the distribution function, = 0, 
corresponds to isotropic initial pressures with Ho = rro = 0. 


Eqs. (23) and (24) reduce to 


n 


e = -- (e + P + n - tt) , 
r 

n . n TT 

— =-onn-H Ahtt — 

Tn r T T 


TT 

T-jr 


4/3^ 
3 r 


1 


TT 2 n 

-1" — 

T 6 T 


(60) 

(61) 

(62) 


where tt = —t^tt'^'^. Note that in this case the term 
involving the vorticity tensor in Eq. (24) vanishes and 
hence has no effect on the dynamics of the fluid. Also 
note that the first terms on the right-hand side of 
Eqs. (61) and (62) corresponds to the first-order terms 
PnO and 2/3^0-'^'^, respectively, whereas the rest are of 
second-order. 

We simultaneously solve Eqs. (60)-(62) assuming an 
initial temperature of Tq = 300 MeV at the initial proper 
time To = 0.5 fm/c, with relaxation times Teq = Tq = 
Ttt = 0.25 fm/c corresponding to initial 77 /s = l/47r. 
We solve these equations for two different initial pres¬ 
sure configurations: ^0 = 0 corresponding to an isotropic 
pressure configuration ttq = Bq = 0 and ^0 = 100 cor¬ 
responding to a highly oblate anisotropic pressure con¬ 
figuration. The anisotropy parameter ^ is related to the 
average longitudinal and transverse momentum in the 
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T [fm/c] 


FIG. 3: (Color online) Same as Fig. 2 except here we consider 
Bose-Einstein statistics. 



r [fm/c] 


FIG. 4: (Color online) Same as Fig. 2 except here we consider 
Fermi-Dirac statistics. 

local rest frame as: ^ = \{Pt)/{Pl) ~ where (• • •) = 
/ d^pfo{\/px + (1 + OPl) A is a temperature-like 

scale which can be identified with the temperature of 
the system in the isotropic equilibrium limit. For par¬ 
ticle mass, we consider m = 300 MeV which roughly 
corresponds to the constituent quark mass. We solve 
Eqs. (60)-(62) with transport coefficients obtained using 
the Grad’s 14-moment method, Eqs. (25)-(29), as well as 
the Chapman-Enskog method, Eqs. (32)-(36). 

In Figs. 2 - 7 we show the time evolution of the bulk 
viscous pressure times proper-time (top) and the pres¬ 
sure anisotropy VlIVt = (P -I- 11 — 7r)/(P -|- 11 -|- 7r/2) 
scaled by that obtained using exact solution of the RTA 
Boltzmann equation (bottom) for three different calcula¬ 
tions: the exact solution of the RTA Boltzmann equation 



T [fm/c] 

FIG. 5; (Color online) Same as Fig. 2 except here we take 
(o ~ 100 corresponding to ttq = 2.86 GeV/fm^ and Ffo = 
0.138 GeV/fmA 



T [fm/c] 

FIG. 6: (Color online) Same as Fig. 5 except here we consider 
Bose-Einstein statistics. 


[33] (red solid line), second-order viscous hydrodynamics 
obtained using the Chapman-Enskog method, Eqs. (32)- 
(36), (blue dashed line) and using the Grad’s 14-moment 
approximation, Eqs. (25)-(29), (green dotted line). Fig¬ 
ures 2 and 5 are for Boltzmann statistics. Figs. 3 and 
6 are for Bose-Einstein statistics, and Figs. 4 and 7 are 
for Fermi-Dirac statistics. Figs. 2-4 correspond to an 
isotropic initial condition (^o = 0), while Figs. 5-7 cor¬ 
respond to a highly oblate anisotropic initial condition 
(^0 = 100 ). 

From Figs. 2 - 7, we observe that, compared to 
Grad’s 14-moment approximation, the transport coeffi¬ 
cients obtained using the Ghapman-Enskog method does 
a marginally better job in reproducing the Pl/Pt ob- 




































T [fm/c] 

FIG. 7: (Color online) Same as Fig. 5 except here we consider 
Fermi-Dirac statistics. 

tained using the exact solution of the RTA Boltzmann 
equation. On the other hand, the result for rll ob¬ 
tained using the Chapman-Enskog method shows better 
agreement with the exact solution of the RTA Boltzmann 
equation than the Grad’s 14-moment method. 

VII. CONCLUSIONS AND OUTLOOK 

In this paper we expressed the transport coefficients 
appearing in the second-order viscous hydrodynamical 
evolution of a massive gas using Bose-Einstein, Boltz¬ 
mann and Fermi-Dirac statistics for the equilibrium dis¬ 
tribution function and Grad’s 14-moment approximation 
as well as the method of Chapman-Enskog expansion for 
the non-equilibrium part. The second-order viscous evo¬ 
lution equations are obtained by coarse graining the rel¬ 
ativistic Boltzmann equation in the relaxation-time ap¬ 
proximation. We also obtained the ratio of the coeffi¬ 
cient of bulk viscosity to that of shear viscosity, in terms 
of the speed of sound, for classical and quantum statis¬ 
tics as well as for the QGP. We then considered the spe¬ 
cific case of a transversally homogeneous and longitudi¬ 


nally boost-invariant system for which it is possible to 
exactly solve the RTA Boltzmann equation [33]. Using 
this solution as a benchmark, we compared the pressure 
anisotropy and bulk viscous pressure evolution obtained 
by employing both the Chapman-Enskog method as well 
as the Grad’s 14-moment method. We demonstrated that 
the Chapman-Enskog method is in better agreement with 
the exact solution of the RTA Boltzmann equation com¬ 
pared to the Grad’s 14-moment method. We found that, 
while both methods give similar results for the pressure 
anisotropy, the Chapman-Enskog method better repro¬ 
duces the exact solution for the bulk viscous pressure 
evolution. 

At this juncture, we would like to clarify that we have 
used the exact solution of the Boltzmann equation, in the 
relaxation-time approximation, as a benchmark in order 
to compare different hydrodynamic formulations. The 
relaxation-time approximation is based on the assump¬ 
tion that the collisions tend to restore the phase-space 
distribution function to its equilibrium value exponen¬ 
tially. Although the microscopic interactions of the con¬ 
stituent particles are not captured in this approximation, 
it is reasonably accurate to describe a system which is 
close to local thermodynamic equilibrium [55]. Looking 
forward, it will be interesting to determine the impact of 
the quantum transport coefficients, obtained herein, in 
higher dimensional simulations. Moreover, it would also 
be instructive to see if the second-order results derived 
herein could be extended to obtain third order transport 
coefficients for quantum statistics [17]. We leave these 
questions for a future work. 
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